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ABSTRACT 

Precise and accurate measurements of neutron star masses and radii would provide 
valuable information about the still uncertain properties of cold matter at supranuclear 
densities. One promising approach to making such measurements involves analysis of 
the X-ray flux oscillations often seen during thermonuclear (type 1) X-ray bursts. 
These oscillations are almost certainly produced by emission from hotter regions on 
the stellar surface modulated by the rotation of the star. One consequence of the rota- 
tion is th at the oscilla tion should appear earlier at higher photon energies than at lower 
energies. iFordl ()l999l ) found compelling evidence for such a hard lead in the tail oscil- 
lations of one type 1 burst from Aql X-l. Subsequently, M uno. Ozel fc Chakrabartvl 
( 2003]) analyzed oscillations in the tails of type 1 bursts observed using RXTE. They 
found significant evidence for variation of the oscillation phase with energy in 13 of the 
51 oscillation trains they analyzed and an apparent linear trend of the phase with en- 
ergy in six of nine average oscillation profiles produced by folding the energy-resolved 
oscillation waveforms from five stars and then averaging them in groups. In four of 
these nine averaged energy-resolved profiles, the oscillation appeared to arrive earlier 
at lower energies than at higher energies. Such a trend is inconsistent with a simple ro- 
tating hot spot model of the burst oscillations and, if confirmed, would mean that this 
model cannot be used to constrain the masses and radii of these stars and would raise 
questions about its applicability to other stars. We have therefore re-analyzed individ- 
ually the oscillations observed in the tails of the four type 1 bursts from 4U 1636—536 
th at, when ave r aged, provided the strongest evidence for a soft lead in the analysis 
bv lMuno et al.1 (|2003lh We have also analyzed the oscillation observed during the su- 
perburst from this star. We find that the data from these five bursts, treated both 
individually and jointly, are fully consistent with a rotating hot spot model. Unfortu- 
nately, the uncertainties in these data are too large to provide interesting constraints 
on the mass and radius of this star. 

Key words: equation of state - relativistic processes - stars: neutron - X-rays: 
binaries - X-rays: burst 



1 INTRODUCTION 

The central regions of neutron stars contain highly de- 
generate matter with densities up to several times nu- 
clear saturation density. The properties of this matter 
cannot be studied in terrestrial laboratories but can be 
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probed by determining the masses and radii o f neu- 
tron stars (see, e.g., iLattimer fc Prakashl 120071 ). One 
method proposed to obtain this information involves fit- 
ting detailed models to the waveforms of the oscilla- 
tions often seen during type 1 ( thermonuclear) X-ray 
bursts (see, e.g., IStrohmaverl 1 19921: iMiller fc Lamb! 1 1998 ; 
Braie et all l200d ; IWeinberg et al.1 l200ll: iNath et all 12003; 
Strohm avcr 200l; ICadeau et al. 2005; Bhattach arvva et al.1 
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l2005l ; ICadeau et al.ll2007l ; lBogdanovll2013l ; see lWattsll2012l for 
a recent overview of the properties of these bursts). These 
oscillations are thought to be produced by rotational modu- 
lation of the emission from a h ot spot that rotates at o r close 
to the stellar spin frequency (Str ohmaver et all ^996). The 
amplitude of the oscillation and its harmonic content are 
affected by special relativistic Doppler boosts and aberra- 
tion and gravitational light deflection, which depend on the 
mass and radius of the star. The mass and radius of the star 
can therefore be determined by fitting sufficiently accurate 
waveform models to waveform data of sufficient quality. 

A strong prediction of the rotating hot spot model is 
that the oscillation should appear earlier at higher photon 
energies than at lower energies. This is most easily seen in an 
extreme case. The spectral shape of the emission from burst 
atmospheres typically has a shape close to the shape of a 
Bose-Einstein spectrum (see, e.g., iBoutloukos et alTl201(t ). 
At energies well above the peak of the spectrum, the ob- 
served maximum of the burst oscillation will occur at the 
stellar rotational phase that maximizes the blueshift of the 
hot spot relative to the observer, which happens when the 
spot is near the approaching limb of the star. In contrast, at 
energies near the peak of the spectrum the maximum of the 
oscillation will occur later, w hen the pro jected area of the 
spot is close to its maximum. iFordl (|l999T l found just such a 
hard lead in Rossi X-ray Timing Explorer (RXTE) obser- 
vations of the oscillation during the tail of a type 1 burst 
fro m Aql X-l. 

iMuno. Ozel fc Chakrabartvl l|2003l ) subsequently ana- 
lyzed 51 oscillation trains in RXTE data from the tails of 
type 1 bursts from six neutron stars. They folded each of 
these oscillation trains to produce oscillation profiles and 
then measured the phases of these profiles using linear least- 
squares fits of sinusoids to the oscillation profiles. Although 
the uncertainties of the oscillation phases were relatively 
large, they found significant evidence for phase variations 
with photon energy in 13 of the 51 oscillation trains. Seven 
of these 13 trains showed random phase variations while the 
remaining six showed either no significant lags or a roughly 
line ar phase trend corr esponding to a hard lag. 

IMuno et all ((2003) then grouped all the energy-resolved 
oscillation profiles from a given star for each RXTE PCA 
gain epoch and averaged the profiles within each of the re- 
sulting nine groups. They found that five of the nine aver- 
aged energy-resolved oscillation profiles showed significant 
phase variations with energy. The phase variations in one of 
these five profiles appeared to be random while the phase 
variations in the other four were consistent with a linear 
increase in phase with energy. The apparent hard lags in 
these four profiles are inconsistent with a simp le rotating 
hot spot model of burst oscillations. IMuno et al.l speculated 
that these hard lags might be caused by Comptonization. 
If confirmed, statistically significant hard lags would mean 
that a simple rotating spot model of burst oscillations can- 
not be used to constrain the masses and radii of these stars 
and would raise questions about the applicability of this 
model to other stars. 

Here we re-analyze individually the oscillations ob- 
served in the tails of four type 1 bursts from 4 U 1636 —536. 
The average oscillation profiles IMuno et al.l {2003) con- 
structed by summing the individual profiles from these 
bursts provided the strongest evidence for a soft lead in their 



burst oscillation analysis. We also analyze the oscillation ob- 
served durin g the superburst from this star. In contrast to 
IMuno et all we find that the data from these bursts, treated 
both individually and jointly, are fully consistent with a ro- 
tating hot spot model. The uncertainties in the RXTE data 
are too large to provide interesting constraints on the mass 
and radius of this star. 

The number of counts in these burst oscillation pro- 
files are so small that the consistency of these profiles with 
the rotating spot model cannot be taken as confirmation 
of this model. However, this consistency does mean that 
these data do not challenge the rotating spot model and 
that fitting this model to burst oscillation data remains a 
viable way to constrain neutron star masses and radii using 
higher quality data from future large-area timing missions 
such as t he Neutron Star Inte rior Composition ExploreR 
(NICER; iGendreau et al. 20121') . the Large Observatory for 



X-ray Timing f£OFr; lMignani et al.ll2012l ; lDel Monte et al.1 
20r j), and the Advanc e d X-ray Timing Array (AXTAR; 



Chakrabartv et al-llioosl ; iRav et aljfeoilfi . 

In Section [2] we describe the X-ray data we use in this 
study and the results of our timing analysis. In Section [3] 
we discuss our results and conclusions. 



2 DATA, ANALYSIS, AND RESULTS 

Our analysis is based on Event Mode data from the RXTE 
Proportional Counter Array (PCA). This mode has a time 
resolution of 1/8192 seconds and is therefore well-suited to 
studying burst oscillations with frequencies of hundreds of 
Hertz. We focused on the 4U 1636—536 type 1 bursts ob- 
served during PCA gain Epoch 4 (the boundaries of the PCA 
energy channels are slightly different in different epochs), 
because the averaged energy-resolved folded profile of the 
oscillations during these bursts that was constructed by 
IMuno et al.l (|2003T ) appeared to show a systematic soft lead, 
contrary to what is predicted by a simple rotating hot spot 
model (see the top middle panel of their Figure 2). We also 
analyzed data from RXTE obse rvations of the hours-long 
super burst from 4U 1636—536 (|Strohmaver fc Markwardg 
I2002T I. which occurred during PCA gain Epoch 5, to de- 
termine whether we could obtain an acceptable joint fit to 
the data on all these bursts. 

There were four type 1 bursts rom 4U 1636—536 during 
PCA gain Epoch 4 that had detectable oscillations; their 
ObsIDs are listed in the first four rows of Table 1. Epoch 4 
ObsIDs 40028-01-08-00 and 40031-01-01-06 each had one ad- 
ditional oscillation train some seconds after the segments we 
analyzed, but these trains were too weak to satisfy our sig- 
nificance criterion and hence we have not included them in 
our analysis. 

We followed IMuno et all (|2003l ) by using five energy 
bins spanning photon energies from 2 and 23 keV, which 
corresponds roughly to PCA energy channels 5 through 53 
(of the 256 total channels). We also followed IMuno et al.l 
by analyzing the oscillations during the tails of the bursts 
rather than during their rising portions, where the changing 
frequencies and rapidly changing amplitudes of the oscilla- 
tion trains could introduce additional complications. 

For each burst, we first searched for the starting time 
and duration of the segment of the burst that maximized 
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Table 1. Properties of the 4U 1636—536 burst segments analyzed 



ObsID 


Start time 


Duration (s) 


Mean Frequency (Hz) 


40028-01-06-00 


171611735.014 


3.260 


580.3985 


40028-01-08-00 


172366987.063 


1.301 


580.5576 


40030-03-04-00 


172431061.302 


1.926 


580.4122 


40031-01-01-06 


172609556.020 


2.645 


581.1129 


50030-02-08-01 


225479593.957 


149.251 


581.9692 



the significance o f the oscillating s ignal. Searches using 
the Leahy power (|Leahv et a l. 1983) or, equivalently, the 
Z\ st atistic (|Buccheri et a l. 1983; S trohmaver fc M arkwardt 
Il999l ). and the deviation of \ 2 f° r the best folded profile 
from its value assuming no oscillation all yielded equivalent 
results. For the four type 1 bursts, we explored durations 
ranging from 0.5 to 6 seconds in steps of 1/16 of a second, 
with starting times 1/16 of a second apart throughout the 
burst tail. For the superburst, we tried durations between 
50 and 250 seconds in steps of 4 seconds throughout the 
portion of the burst where X-ray flux oscillations were ob- 
served. For the first type 1 burst and the superburst, we also 
tried fits with a frequency that varies linearly in time; the 
results were indistinguishable from those for our fits assum- 
ing a constant fre quency. We then used epoch folding (see, 
e.g., lDavieslll990r i to determine the frequency of the oscil- 
lation during each burst segment and to construct a folded 
oscillation profile for each set of energy channels. The start 
time, duration, and oscillation frequency for each of the five 
burst segments we analyzed are listed in Table [T] 

We determined the zero phase (maximum flux) of the 
bolometric oscillation profile and the oscillation profile in 
each set of energy channels by fitting a sinusoidal model to 
each of these profiles. The parameters in this model were the 
phase, amplitude, and DC level. We estimated the la un- 
certainties in these parameters using the standard approach. 
Namely, we began with the overall best-fit values of the three 
parameters. We then varied one parameter from its best-fit 
value while minimizing \ 2 with respect to the other two pa- 
rameters, until we found the value of the first parameter that 
gave a A\ 2 of 1 relative to the x 2 f° r the best-fit parameter 
values. We used this value of the first parameter to define 
its lcr uncertainty region. We then repeated this procedure 
for the other two parameters. 

Following iMuno et al.l l|2003i ). we used the parameter 
values in our fitted sinusoidal models to determine the phase 
of the oscillation profile in each set of energy channels rel- 
ative to the phase of the bolometric oscillation profile. We 
also computed the relative phases by cross-correlating the 
profiles. The results were consistent with the results from 
the sinusoidal fits. Both methods confirmed the hard lead 
found bv lFordl (|l999l ) for the oscillation in a segment of the 
tail of a burst from Aql X-l. We also studied earlier por- 
tions of that burst but found large data gaps that render 
any conclusions concerning these portions untrustworthy. 

Figure [T] shows our results for the phases of the oscil- 
lation profiles as a function of photon energy. These results 
are insensitive to small changes in the start times and dura- 
tions of the data segment and the oscillation frequencies. In 
this figure, a negative lag implies that the maximum of the 
oscillation in that energy range arrived, on average, before 



the maximum of the oscillation in the bolometric profile. A 
positive lag means the opposite. As explained earlier, a sim- 
ple rotating spot model predicts an increasingly negative lag 
with increasing energy. 

Figure [T] also shows the results of a fit to these data 
of the model waveform produced by a uniform circular hot 
spot fixed to the rotating star, computed usi ng the algo- 
rithms described briefly in lLamb et al.l (|2009al rbh and more 
extensively in Lo et al. (2013, in preparation). In fitting this 
model, the mass, radius, observer latitude, and distance to 
the star were kept fixed for all five bursts, but the latitude of 
the spot center, the spot angular radius, and the spot color 
temperature were allowed to be different for each burst. In 
order to extract maximal information from the data, we used 
standard Bayesian techniques to compare our model predic- 
tions with the data, which were broken into 16 phases and 
25 energy channel groupings per phase for each burst. The 
fit shown is for a star with a mass of 1.5 Mq and a ra- 
dius of 8.8 km, an observer inclination relative to the stellar 
spin axis of 56°, and a distance of 6.5 kpc (consistent with 
the distance of (6.0 ± 0.5)(M/1.4 M©) 1/2 kpc estimated by 
iGallowav et"ail 2006). A wide range of other stellar masses 
and radii give comparably good fits to these data, so we are 
unable to derive interesting constraints on the properties of 
this neutron star using these data. 

A simple rotating hot spot model is fully consistent with 
the 4U 1636—536 data that we analyzed. The fits of this 
model go through 76% (1 9 of 25) of the one-sigma error re- 
gions. The Cash statistic (|Cashll 19791 ; note that this statistic 
is an analog of \ 2 f° r small numbers of counts and asymp- 
totes to x 2 f° r large numbers of counts) for the joint fit is 
1859 for 1850 degrees of freedom. If we treat this as a \ 2 
distribution, then if the model is correct we will find \ 2 °f 
1859 or larger for 1850 degrees of freedom in approximately 
44% of realizations. This is therefore an excellent fit, espe- 
cially because the data shown for each burst were obtained 
by folding the oscillation waveform over several seconds of a 
burst that itself lasted only a few seconds, a procedure that 
undoubtedly smeared the oscillation profile. The uncertain- 
ties in the lag estimates are sufficiently large that occasional 
apparent soft leads (such as the one in the leftmost panel) 
are consistent with being statistical fluctuations. The su- 
perburst, which has far more counts than any of the type 1 
bursts, does appear to show a hard lead, in accordance with 
the prediction of a simple rotating spot model. 
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Figure 1. Phase lags for oscillations in the tails of five X-ray bursts from 4U 1636—536 (data points) and a joint fit of a simple rotating 
circular hot spot model to all the data (broken lines), showing that this model is fully consistent with the data. The phases shown for 
the data points and the model points are both relative to the zero phase of the observed bolometric oscillation profile, which is defined 
by the maximum flux of the profile. Each panel shows the data and ObsID for one of the five bursts. The horizontal error bars indicate 
the approximate energy range covered by each channel while the vertical error bars show the la uncertainties determined by fitting a 
sinusoid to the data. In fitting the rotating spot model to these data, the mass, radius, observer latitude, and distance to the star were 
kept fixed for all five bursts, but the spot latitude, spot radius, and spot color temperature were allowed to vary from burst to burst. See 
the text for further details. 



3 DISCUSSION 

We have individually analyzed the X-ray flux oscillations 
during four type 1 X-ray bursts from 4U 1636—536 observed 
using RXTE during its PCA gain Epoch 4 and the super- 
burst from this star. We find that the variations of the phases 
of these oscillations with photon energy are fully consistent 
with a rotating hot spot model, whether the data from each 
burst is fit individually or jointly. Unfortunately, the uncer- 
tainties in these data are too large to provide interesting 
constraints on the mass and radius of this star. 

Our results for the variations of the phases of the os- 
cillation pr ofiles with photon energy differ from the results 
reported bv lMuno et alj |2003l ). which appeared to be incon- 
sis tent with a s imple rotating hot spot model. The analysis 
bv lMuno et al.l diffe red from ours in several respects. In par- 
ticular, Muno et alj averaged the energy-resolved profiles of 
the oscillations observed during the type 1 X-ray bursts from 
4U 1636-536 during RXTE PCA gain Epoch 4 and then an- 
alyzed them, whereas we analyzed the energy-resolved pro- 
files of the oscillations observed during each burst individu- 
ally and then jointly. 

Averaging the oscillation profiles before analyzing them 
is not the optimal procedure, because the variations of the 
phases of the oscillation profiles with energy differ from burst 



to burst in this data set, as Figure [T] shows. This is to be 
expected, because the properties of the hot spot (such as its 
radius or inclination fro m the spin axis) can change from 
burst to burst. Indeed, iMuno et al.l found that the varia- 
tions with energy of the phases of the oscillation profiles in 
the 4U 1636—536 bursts changed with time (the dependence 
during PCA gain Epochs 3, 4, and 5 differed from one an- 
other). This suggests that the oscillation profile does vary 
from burst to burst. 

Our most important result is that the X-ray flux oscilla- 
tions during the four type 1 X-ray bursts from 4U 1636—536 
that we analyzed and the superburst from this star are fully 
consistent with a simple rotating hot spot model, whether 
the data from each burst are fit individually or jointly. Thus, 
these data do not call into question the use of such a model 
to fit X-ray burst oscillation data and thereby constrain the 
masses and radii of individual neutron stars. At the same 
time, it is clear that more and better data will be required 
to obtain tight constraints. 

Proposed space missions that will obtain much bet- 
ter timing data include NICER, which will focus on 
deep observation s of X-ra y emitting millisecond pulsars 
i|Gendreau et al.l |2012| ; see iBogdanovl 2013 for constraints 
recently obtained using XMM-Newton data); LOFT, which 
will have a collecting area more than an order of mag- 
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nitude larger than the R XTE PCA (see iMignani et al.l 
l2012l ; lDel Monte et alj|2012[ ); and AXTAR, which would also 
have a collecting area much larger than the RXTE PCA 
l|Chakrabartv et alJliooj ; iRav et alJlioilh . Our simulations 
of the precision of the phase lag measurements that could be 
achieved using LOFT show that it could achieve a precision 
of ~0.01 cycles, compared to the ~0.05 precision achieved 
using RXTE. A full Bayesian analysis of the constraints on 
the masses and radii of neutron stars that could be achieved 
by fitting waveform models to LOFT observations of burst 
oscillations (Lo et al. 2013, in preparation) indicates that 
tight constraints can be achieved for systems that have a fa- 
vorable geometry. Thus, there is a good prospect that these 
missions will provide neutron star mass and radius estimates 
precise enough to tightly constrain the properties of cold 
supranuclear matter. 
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